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ABSTRACT 


The  effect  of  gravity  forces  on  turbulence  development 
is  investigated  both  experimentally  and  numerically. 

In  the  experimental  investigation,  helium  gas  is 
released  into  an  air  stream.  The  experiment  is  carried  out  in 
a  low  speed,  low  turbulence  wind  tunnel.  The  He  is  released 
through  an  area  source.  Two  configurations  are  investigated  : 

He  released  from  the  top  of  the  section,  and  He  released  from 
the  bottom  of  the  section.  Three  different  free  stream  turbu¬ 
lence  levels  are  investigated  at  two  free  stream  velocities. 

Numerically,  the  conservation  of  mass  equation,  written 
for  a  system  of  two  species,  is  solved.  The  turbulent  diffusion 
terms  are  modelled  with  a  gradient  flux  model,  constant  coef¬ 
ficient  of  diffusivity.  Comparisons  are  made  with  the  experi¬ 
mental  data,  comparing  effects  of  changing  turbulent  intensity, 
free  stream  velocity  and  buoyancy. 


LIST  OF  SYMBOLS 


mass  fraction  of  gas 
frequency  of  laser  light 

stretching  function  in  x,  acceleration  due  to  gravity 

stretching  function  in  y 

current  through  hot  wire 

length  of  hot  wire 

chemical  production 

time 

velocity  vector 
molecular  diffusivity 
wavenumber  vector,  diffusivity 
potential  energy 
resistance  through  hot  wire 
temperature 

ratio  of  density  of  air  to  density  of  gas 

boundary  layer  thickness,  function  of  concentration 
i ncrement 

dummy  variable 

wavelength  of  laser  light,  weight  for  upwinding 
total  mass  concentration 
mass  concentration  of  gas 

Superscripts 

fluctuating  quantities,  derivative 

mean  quantities,  transformed  coordinate  system 
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1.  INTRODUCTION 

The  study  of  the  dispersal  of  pollutants  is  of  major 
concern  to  all  factions  of  society  :  governments,  industries, 
communities,  and  to  each  individual  personally.  With  the  ex¬ 
panded  use  of  petroleum  products,  including  LPG  and  large  chemi¬ 
cal  production  plants,  interest  has  recently  been  escalated 
regarding  the  probable  outcome  of  large  accidents  with  heavy 
gases.  In  the  past,  very  approximate  methods,  usually  loaded 
with  empiricism,  were  used  to  predict  the  average  concentration 
field  of  a  pollutant.  This  is  satisfactory  if  it  is  applied 
to  mildly-toxic  or  non  flammable  emissions,  and  usually  very 
large  safety  factors  are  applied  to  cover  the  weaknesses  in 
the  schemes.  But  in  this  day  and  age,  where  large  amounts  of 
highly  toxic  and  flammable  gases  are  being  stored  and  transported, 
more  accurate  prediction  methods  are  needed.  To  try  to  meet 
this  need,  many  new  models  have  been  proposed,  based  on  the 
physics  of  the  problem.  That  is,  an  attempt  to  solve  the 
governing  equations  of  the  flow  is  made.  Since  turbulence  on 
all  scales  is  present  in  the  flow  field,  assumptions  and  simpli¬ 
fications  must  be  made  to  make  these  equations  soluble.  Thus, 
the  results  of  such  studies  are  no  more  accurate  than  the  empi¬ 
rical  models.  The  models  also  vary  in  the  results,  one  to  the 
other,  as  in  each,  the  simplifications  are  made  in  different 
ways  and  in  varying  degrees.  A  major  problem  with  these 
studies  is  the  scarcity  of  good  experimental  data,  especially 
for  the  dispersion  of  heavy  gases. 

An  extensive  review  of  models  and  methods  used  in 
pollutant  dispersal  for  non-buoyant  gases  can  be  found  in  Slade 
(Ref.  1).  For  a  review  of  models  as  particularly  applied  to 
heavy  gases,  reference  can  be  made  to  the  paper  by  Raj 
(Ref.  2). 

In  the  present  study,  a  simple  flow  field  is  studied 
experimentally  to  provide  data  to  increase  the  understanding  of 
the  phenomena  of  heavy  gas  dispersal.  These  data  are  then  com¬ 
pared  to  a  simple  numerical  scheme,  using  a  constant  coefficient 
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of  turbulent  diffusivity,  with  the  purpose  of  determining  the 
validity  of  such  a  model. 

The  experimental  investigation  consists  of  the  release 
of  He  gas  from  an  area  source  in  a  wind  tunnel.  Both  negatively 
and  positively  buoyant  gas  releases  are  studied  by  releasing  the 
gas  from  either  the  top  of  the  section  or  the  bottom  of  the 
section.  Two  average  flow  velocities  at  three  turbulent  inten¬ 
sities  are  studied.  The  release  is  such  that  the  flow  field  can 
be  considered  two  dimensional. 

The  numerical  study  consists  of  the  solution  of  the 
diffusion  equation  by  approximate  factorization  after  the  manner 
of  Beam  and  Warming  (Ref.  3).  The  velocity  is  assumed  to  be 
unaffected  by  the  gas,  eliminating  the  need  to  solve  the  momen¬ 
tum  equation.  As  indicated  in  the  previous  paragraph,  the  field 
is  considered  two  dimensional.  Diffusion  is  non-negl igible  in 
both  spatial  directions,  but  convection  in  the  vertical  direc¬ 
tion  (perpendicular  to  the  direction  of  the  flow)  is  neglected. 


2.  THEORY 


In  the  case  of  the  dispersal  of  a  pollutant  into  the 
atmosphere,  one  must  consider  a  flow  field  with  two  gaseous 
species.  The  equation  describing  the  diffusion  of  one  species 
into  another  can  be  written  (Ref.  4,  p  557)  : 


3p  a 
— -  + 

at 


$*paV  = 


V*p 


dab7c 


+  r, 


where 


(2.1) 


r^  =  chemical  production 

=  mass  concentration  of  gas 
p  =  total  mass  concentration 

c  *  mass  fraction  of  gas 

Dab  =  molecular  diffusivity  of  gas  into  air. 

In  the  case  under  study,  there  are  no  chemical  reactions,  so 

rA  =  PA’  the  mass  concentration,  can  be  expressed  as  cp^, 

where  p  is  the  density  of  the  gas  and  c  is  the  mass  fraction 

of  gas.  If  p_  is  the  density  of  air,  p  =  c(p  -p  1  +  p  .  Then 
o  g  a  a 

equation  (2.1)  can  be  written  as  : 

dab  —  (2-2> 

AB 

(Note  that  the  usual  summation  convention  is  inferred  by  repeated 
i ndi ces ) . 


If  it  is  further  assumed  that  : 
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c  =  c  +  c ' 


P  -  Y(Pg-Pa)  +  Pa 


pa/og 
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and  equation  (2.2)  is  averaged  in  time,  we  obtain  : 


i£  t  _L.  cili+  -i-  cvT  =  (1-0.)  -5-  D  !£- 

^  w  * 


3 1  3xi 


3Xi 


.  AB  . 
3x^  3  X .. 


(2.3) 


Equation  (2.3)  can  be  solved  for  the  composition  ci 
the  mixture  if  the  velocity  field  and  the  turbulent  dispersion 
terms  are  known.  Normally,  the  equation  of  motion  is  solved 
for  the  velocity  field  and  the  turbulent  dispersion  terms  are 
modelled  in  some  way.  In  this  work,  it  is  assumed  that  the  velo 
city  is  not  changed  by  the  presence  of  the  gas,  so  it  can  be  a 
known  input  to  equation  (2.3),  and  it  is  not  necessary  to  solve 
the  equation  of  motion.  This  simplifies  the  problem  immensely, 
as  then  it  is  necessary  only  to  solve  one  equation,  with  a 
model  for  the  turbulent  dispersion  terms. 


EXPERIMENTAL  INVESTIGATION 


3 . 1  Experimental  method 

3.1.1  VeJ_oc  i  ty_meas^remejnts_ 

A  complete  explanation  of  the  LDV  techniques  and 
associated  bibliography  can  be  found  in  reference  5. 

Here  the  basic  principles  of  this  measurement  technique  and 
the  motivation  of  its  employment  in  turbulent  measurements 
will  be  briefly  recalled. 

The  laser  doppler  method  is  based  on  the  detection 
of  the  doppler  frequency  of  laser  light  scattered  from  small 
solid  or  liquid  particles  moving  with  the  main  flow. 

The  relation  between  the  velocity  of  the  particles 
illuminated  by  the  laser  beam,  of  frequency  f  and  wavelength 
A,  and  the  doppler  frequency  fQ  of  the  scattered  light  is 
linear.  That  is, 

Vfs  “  U/A(KS-K.)  (3.1) 

where 

f  =  frequency 
K  =  wave  number  vector 
U  =  velocity  vector 
subscript  i  =  incident  light 
subscript  s  =  scattered  light 

With  the  LDV  it  is  also  possible  to  measure  simulta¬ 
neously  two  velocity  components  in  a  turbulent  flow.  In  fact, 
if  two  pairs  of  orthogonal  laser  light  beams  are  focused  on  the 
same  point,  a  spatial  system  of  internally  orthogonal  interfe¬ 
rence  fringes  is  formed.  If  one  pair  of  beams  is  shifted  in 
frequency,  for  example  by  passage  through  a  light  frequency 
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shift  unit  such  as  a  Bragg  cell,  the  probe  volume  is  composed 
of  the  superposition  of  a  stationary  fringe  system  and  a  moving 
one,  which  has  a  fringe  velocity  proportional  to  the  frequency 
difference  of  the  two  shifted  beams.  A  particle  passing  through 
this  volume  scatters  light  of  two  different  frequencies  corres¬ 
ponding  to  the  components  of  the  particle  velocity  with  respect 
to  fringe  planes.  Therefore,  it  is  possible  to  separate  the 
frequency  corresponding  to  each  velocity  component  by  filtering 
the  signal  with  two  different  filters.  The  necessary  conditions 
to  achieve  this  result  are  : 

1.  no  overlapping  in  the  frequencies  produced  by  the  two 
fringe  systems,  and 

2.  a  detectable  frequency  difference. 

From  these  considerations  it  is  clear  that  a  twin  LDV 
with  Bragg  cell  in  one  arm  is  a  very  suitable  system  for  two 
dimensional  velocity  measurements  in  a  turbulent  flow  because 
of  its  non-i ntrusi vi ty ,  which  is  a  very  stringent  requirement 
for  the  measurements  of  turbulent  quantities,  and  because  the 
two  measurements  can  be  taken  simultaneously. 

3.1.2  £once^ntr£ti^on  measurements 

Helium  concentration  is  measured  with  a  probe  of  the 
kind  introduced  and  tested  by  Olivari  &  Colin  (Ref.  6).  The  work¬ 
ing  principle  is  the  following.  If  a  metallic  wire,  heated  by  an 
electrical  current,  is  placed  in  a  flow,  the  voltage  difference 
between  its  extremes  is  sensitive  to  the  variation  of  velocity 
and  heat  transfer  properties  of  the  flow  itself,  and  thus  to 
its  composition,  if  the  medium  is  a  mixture  of  two  gases.  There¬ 
fore,  in  principle,  it  is  sufficient  to  separate  these  indepen¬ 
dent  variables  to  measure  the  variation  of  one  of  them,  e.g. 
in  order  to  determine  the  concentration,  one  can  avoid  the  velo¬ 
city  dependence  by  placing  the  wire  in  a  region  of  uniform 
velocity.  To  create  the  uniform  velocity  region,  the  probe 
has  a  constant  section  followed  by  an  expansion.  The  wire  is 


placed  in  the  center  of  the  expansion  region,  where  the  flow 
behaves  like  the  potential  core  of  a  jet.  The  mass  flow  in  the 
probe  is  controlled  by  a  sonic  hole  driven  by  a  vacuum  pump. 

In  such  a  manner  the  mass  flow  is  a  function  only  of  the  con¬ 
centration  (and  eventually  temperature)  variation  since  this 
parameter  affects  the  speed  of  sound. 

If  the  hot  wire  is  placed  in  a  region  of  relatively 
high  velocity,  which  minimizes  its  sensitivity  to  Re  variation, 
the  simplified  equation  for  the  hot  wire  heat  transfer  is  : 

i 2R/A  =  lk (T  -T  .  )  *  AlAAt  (3.2) 

'  w  mix'  v  ' 

and  if  the  hot  wire  is  inserted  in  a  resistance  bridge,  the 

voltage  across  the  bridge  is  : 

v  =  Rb(1+1‘)  *  RBi  (3.3) 

Thus 

1/2 

V  =  Rb(LXAAT/R)  (3.4) 

From  this  relation  it  is  evident  that  there  exists  a 
dependence  of  the  voltage  wire  output  on  both  thermal  conduc¬ 
tivity  and  temperature,  so  it  is  again  necessary  to  carry  out 
a  variable  separation,  that  is,  running  at  a  constant  flow 
temperature.  Also  evident  is  the  advantage  of  the  employment 
of  a  hot  wire  with  very  low  electrical  resistance  in  order  to 
increase  the  sensitivity. 

3 . 2  Test  facilities  and  instrumentation 

The  measurements  are  carried  out  in  a  low  speed  wind 
tunnel  (VKI  L-3)  whose  features  are  :  blower-driven,  open-return 
circuit,  25/1  convergence  contraction  ratio,  temperature 
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controlled,  free  stream  turbulence  level  of  0.1%,  cross  section 
of  200  mm  *  200  mm  and  a  test  section  2  m  long. 

A  new  test  section  was  designed  in  order  to  generate 
a  uniform  horizontal  flow  with  turbulence  level  varying  in  the 
range  2-10%  and  a  density  stratification  in  the  vertical 
di recti  on . 


3.2.1  FI owsys tern 

3. 2. 1.1  Turbul enceproducti on 

The  different  levels  of  turbulence  are  obtained  by 
placing  different  kinds  of  turbulence  generators  and  manipula¬ 
tors  near  the  entrance  of  the  test  section.  The  choice  of  these 
turbulence  manipulators  is  based  upon  the  results  and  informa¬ 
tion  contained  in  the  work  reported  in  reference  7,  dealing 
with  the  interaction  of  free  stream  turbulence  with  screens  and 
grids. 


Table  3.1  specifies  the  manipulators  selected,  three 
different  screens  and  one  perforated  plate,  and  summarizes 
some  of  their  relevant  characteristics.  The  table  also  gives 
the  pressure  drop  coefficient  k,  the  mesh  Reynolds  number,  Rem 
and  the  Reynolds  number,  Re^,  based  on  the  wire  diameter  of 
the  screens,  for  a  free  stream  velocity  of  3  m/s. 

The  solidity,  o,  is  always  less  than  45%  to  avoid  the 
emergence  of  large  scale  instability  associated  with  high  soli¬ 
dity  devices.  The  Red  is  always  more  than  30,  which  is  the 
critical  range  in  which  small  changes  in  U  or  upstream  distur¬ 
bances  result  in  large  changes  in  the  character  of  the  flow 
downstream  of  the  manipulators. 

For  the  final  six  flow  conditions,  the  manipulators 
were  chosen  to  have  turbulence  intensities  of  16%  (manipulator  1 
8%  (manipulator  2)  and  3%  (manipulator  3). 


3.2. 1.2  Density  gradi ent_generati on 

The  density  gradient  in  the  vertical  direction  is 
obtained  by  injecting  helium  through  a  porous  wall  section.  The 
test  section  can  be  adapted  such  that  the  porous  wall  is  on  the 
top  or  on  the  bottom  side.  When  He  is  released  from  the  top 
wall,  a  stable  stratification  results,  as  all  the  energy 
generated  by  shear  is  used  in  working  against  the  buoyancy 
forces,  which  therefore  produce  a  loss  of  turbulence  energy  in 
addition  to  viscous  dissipation.  When  He  is  released  from  the 
bottom  wall,  however,  an  unstable  stratification  results,  as 
buoyancy  forces  transfer  energy  to  the  turbulence  at  a  rate 
that  is  independent  of  height,  while  transfer  from  the  mean  flow 
is  proportional  to  the  velocity  gradient  and  buoyancy  becomes 
relatively  larger;  the  ratio  of  turbulent  intensity  to  shear 
stress  increases  and  the  eddy  diffusivity  is  larger  than  that 
in  a  constant  density  flow. 

The  helium  is  injected  through  the  porous  section 
with  a  volume  flow  of  4740  1/ h. 

The  porous  wall  is  made  of  plastic  foam  that  allows  the 
gas  to  perspire  with  an  exit  velocity  practically  equal  to  the 
filter  velocity,  1.9  cm/s.  As  the  ratio  between  the  filter 
velocity  and  the  air  velocity  is  of  the  order  of  1/150  or  1/300 
(u  =  3  m/s  and  6  m/s,  respectively),  it  is  possible  to  neglect 
the  perturbation  due  to  the  lateral  injection  of  flow.  The 
reference  system  is  therefore  taken  as  the  main  air  flow 
entering  the  test  section. 

3.2.2  Ex£e r i mental _a£p a ratus 

3. 2. 2.1  Velocity  measurements 

A  Spectra-Physi cs  model  12A  helium-neon  laser  source 
delivering  a  power  of  15  mW  is  employed.  After  a  first  polari¬ 
zation  the  laser  beam  enters  the  transmitting  part  of  the 


N  > 


-  10  - 

optical  system,  which  is  composed  of  the  following  elements  : 
triple-beam-splitter  module,  which  splits  the  beam  into  three 
separate  beams,  two  beams  of  equal  intensity  at  a  distance  of 
50  mm  apart  and  a  central  beam;  a  second  polarizer  for  the 
central  beam;  dual-beam  splitter  that  splits  the  central  beam 
into  two  equal  intensity  beams  with  a  separation  of  20  mm;  the 
double  Bragg  cell  used  to  shift  the  frequency  of  one  of  the 
two  pairs  of  beams;  the  steering  module  to  center  the  shifted 
beams  in  the  plane  orthogonal  to  the  non  shifted  ones;  and, 
last,  a  focusing  lens  with  focal  length  of  333  mm. 

The  receiving  optics  consist  of  a  photomultiplier 
tube  RCA-97326  focused  by  a  zoom  objective.  The  signal  from 
the  photomul tipi ier  is  then  band-pass  filtered  in  two  DO-781-3 
VKI  filters  in  order  to  discriminate  the  two  velocity  components 
(Ref.  8).  Afterwards  the  signals  enter  two  VKI-DO-78-31  counter 
data  processors  which  give  fluctuating  voltage  outputs  directly 
proportional  to  the  velocity  components. 

The  most  important  parameters  defining  this  optical 
system  are  reviewed  in  table  3.2  :  probe  volume  dimensions; 
number  of  figures  inside  it;  fringes  spacing;  and  pinhole 
diameter.  In  figure  3.1  a  schematic  view  of  the  source  with 
the  electronics  is  shown. 

3. 2. 2. 2  Concentration  measurements 

For  concentration  measurements  the  techniques  described 
in  section  3.1.2  are  adopted.  The  miniaturized  gas  concentra¬ 
tion  probe  is  placed  in  the  flow  field.  A  mechanical  device, 
fixed  on  the  external  part  of  the  test  section  wall,  supports 
it  and  allows  its  displacement  in  the  vertical  direction  (ortho¬ 
gonal  to  the  flow).  The  output  of  the  probe  is  fed  to  the  VKI 
series  85  hot  wire  anemometer  and  then  to  an  oscilloscope  to 
visualize  the  signal . 


This  probe  is  not  an  absolute  concentration  sensor 
and  as  such  it  needs  a  previous  calibration  for  the  gas  employed 
in  the  experiment.  Moreover,  since  it  is  very  sensitive  to  the 
ambient  temperature  the  calibration  is  repeated  each  day  just 
before  starting  the  measurements.  The  sensitivity  near  the 
origin  is  0.15%  He  in  air  in  volume  per  millivolt.  The  response 
time  is  of  the  order  of  1  m/s.  Figure  3.2  is  a  schematic  of 
the  equipment  employed. 

3. 2. 2. 3  System_of  data  acqui si tion_and_ reduction 

Concentration  and  velocity  profiles  are  taken  in 
several  sections  along  and  downstream  of  the  surface  gas  source. 
The  main  quantities  to  be  obtained  in  order  to  understand  the 
behavior  of  the  concentration  and  velocity  fields  are  the  mean 
velocity  and  concentration,  the  RMS  for  the  two  velocity  compo¬ 
nents  and  the  concentrati on ,  and  their  respective  spectra.  To 
obtain  all  of  that  it  is  necessary  to  collect  and  store  large 
amounts  of  data  for  each  measurement  point.  For  this  reason 
the  standard  VKI  medium  speed  system  of  data  acquisition  on  the 
PDP  11/34  is  utilized.  The  scheme  is  shown  in  figure  3.3.  All 
the  signals,  those  from  the  LDV  (two  velocity  components)  and 
those  from  the  HWA  (concentration),  before  arrival  in  the  ADC, 
are  low-pass  filtered  at  1.5  kHz  to  eliminate  noise  (1.5  kHz 
is  a  value  consistent  with  the  response  time  of  the  concentra¬ 
tion  probe  and  with  the  value  of  the  data  rate  of  the  LDV). 

The  data  were  digitized  and  stored  on  the  computer  to  be  pro¬ 
cessed  afterwards  using  VKI  transfer  function  programs  available 
on  the  PDP  11/34  and  on  the  VAX. 

3.2.3  Preci_s_i_on  of_th_e_meas_urernents^ 

Following  the  standard  error  analysis  (Ref.  9)  the 
experimental  uncertainty  is  of  the  order  of  3%  for  the  velocity 
measurements  and  of  the  order  of  2%  He  in  air  in  volume  for  the 
concentration  measurements.  For  the  power  spectra  computation 
with  the  FFT  VKI  program  the  frequency  bandwidth  is  6  Hz. 
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3.3  Presentation  of  results 

3.3.1  Me  an_con£e  it t  r  a  t  i  o  n_me  a  s^i  r  ement  s^ 

The  mean  concentration  profiles,  for  the  16%  turbulence 
case,  at  constant  x  positions,  are  shown  in  figures  3.4  and  3.5 
(all  the  data  for  the  six  cases  can  be  found  in  reference  10). 

It  can  be  clearly  seen  from  these  figures  that  in  the  unstable 
case,  dissipation  is  increased,  which  is  the  expected  result. 

The  same  data,  but  plotted  in  terms  of  constant  y 
positions,  are  shown  in  figures  3.6  and  3.7.  At  y  positions 
very  close  to  the  wall,  a  certain  flattening  is  observed.  This 
can  arise  from  the  averaging  effect  of  the  probe  size  and  also 
from  the  effect  of  the  vertical  velocity  component. 

3.3.2  Me an_ve 1 o c measurements 

For  the  mean  velocity  profiles,  the  most  interesting 
observations  come  from  the  comparison  of  the  three  configura¬ 
tions:  reference,  stable  and  unstable.  The  shape  of  the  profile 
itself  is  not  very  sensitive  to  changes  in  parameters  of  the 
turbulent  layer  or  the  mean  external  velocity. 

As  a  typical  example,  compare  the  velocity  profiles 
taken  at  a  distance  equal  to  52.5  cm  downstream  from  the  leading 
edge  of  the  plate,  for  a  turbulence  level  of  3%.  These  profiles 
are  shown  in  figure  3.8.  It  is  evident  that  stratification  in 
density  influences  the  thickness  of  the  boundary  layer. 

3.3.3  Tujrbijl  en£e_i  nten^i  ties_measurement£ 

The  turbulent  intensity  is  defined  as  the  ratio  of 
the  root-mean-square  of  the  fluctuating  part  to  the  mean  part. 
Figure  3.9  shows  the  turbulent  intensity  profiles  of  the  longi¬ 
tudinal  and  vertical  velocity  components  and  of  the  He  concen¬ 
tration  for  the  same  geometrical  and  flow  conditions  as 


13 


previously  discussed.  In  the  unstable  case  the  turbulent 
layers  are  about  two  times  the  ones  of  the  neutral  condition 
in  the  first  50%  of  the  boundary  layer  thickness  while,  outside 
the  boundary  layer  thickness,  they  have  about  the  same  value. 

In  the  stable  stratification,  near  the  wall  (first 
50%  of  6),  the  turbulent  layers  remain  larger  than  the  ones  of 
the  reference  case,  although  much  smaller  than  the  correspond!' ng 
values  for  the  unstable  stratification. 

3.3.4  Turbu/1  en£e_spectra_and_nucro  scales 

The  power  spectral  density,  in  terms  of  frequency  for 
the  fluctuating  part  of  the  concentration  and  of  the  two  velo¬ 
city  components  are  also  computed  with  the  method  of  the  fast 
Fourier  transform.  Details  and  results  of  these  computations 
can  be  found  in  reference  10.  Some  spectra  for  a  specified 
spatial  position  are  shown  in  figure  3.9.  Both  stable  and  un¬ 
stable  cases  are  shown  for  the  low  turbulence  case.  There  is 
not  a  great  deal  of  difference  between  the  two  curves,  although 
it  can  be  said  that  there  is  a  small  increase  in  power  for  the 
unstable  case,  an  expected  trend.  The  difference  is  not,  however, 
large  enough  to  make  any  conclusive  statements.  The  third  curve 
is  for  the  stable  case  at  a  higher  turbulence  level  (16%).  Here 
the  high  frequencies,  representing  the  background  turbulence, 
have  a  higher  power  level,  of  course.  The  other  qualities  of 
the  curve  remain  basically  unchanged. 

3.4  Analysis  of  possible  errors 
in  concentration  measurements 


3.4.1  Probe  si^e 

The  diameter  of  the  probe,  at  its  widest  part,  is 
9  mm.  Very  little  reliability  can  be  placed  on  the  measurements 
at  z  =  0.4  and  0.1  mm,  as  the  disturbance  to  the  flow  is  on 
this  scale.  The  opening  of  the  probe  itself  is  2  mm  in  diameter, 


meaning  that  any  measurement  is  an  average  over  at  least  this 
2  mm  diameter  space.  Close  to  the  wall,  where  the  gradients 
are  expected  to  be  very  steep,  this  may  result  in  serious 
i naccuracies . 

3.4.2  Calibration 

Although  an  error  analysis  showed  that  the  concentra 
tion  can  be  measured  to  within  2%  accuracy  (see  experimental 
section)  inspection  of  figure  3.10  shows  that  a  linear  calibra 
tion  curve  fits  the  data  to  only  about  7%.  Thus,  the  concen¬ 
tration  measurements  cannot  be  expected  to  be  better  than  7% 
accurate. 
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4.  NUMERICAL  SIMULATION 

4 . 1  Basic  equations  and  simplifications 

The  equation  to  be  solved  numerically  is  (refer  to 
equation  2.3)  : 


3C  3  —  3  i  r  . ,  .  3  _  3 c 

— —  +  -  cu.  +  -  c'u.  =  ( 1-a)  -  D ,  „  - 

3t  3X .  1  3xi  3Xi  3xi 


In  a  two  dimensional  flow  field,  this  becomes 


(2.3) 
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The  coordinate  system  is  set  with  x  being  in  the  direction  of 
the  flow  and  y  is  in  the  vertical  direction. 

In  the  actual  experiment  (see  Chapter  3),  the  exit 
velocity  of  the  gas,  through  the  porous  plate,  is  approximately 
0.02  m/s.  An  approximate  calculation  (see  section  4.6.3) 
showed  a  maximum  velocity  contribution  from  gravity  effects  to 
be  of  the  order  of  0.08  m/s.  This  can  be  assumed  equal  to  the 
maximum  value  of  v  in  the  field  (mean  vertical  velocity 
component).  Thus,  the  third  term  in  equation  (4.1)  may  be 
neglected,  yielding  : 
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Ordinarily,  as  mentioned  previously,  the  momentum 
equation  is  solved  for  the  mean  velocity  in  the  horizontal 
direction,  u,  for  input  into  equation  (4.2).  But  if  it  is 
assumed  that  the  presence  of  the  gas  does  not  affect  the  mean 
velocity  field,  u  can  be  considered  a  known  input,  reducing 
the  problem  to  that  of  solving  one  equation  only,  equation  (4.2) 


and  v ' c ' , 


For  the  turbulent  dispersion  terms,  u'c' 
the  gradient  flux  model  is  used,  i.e.. 


u'c'  *  -K  —  (4.3) 
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where  K  and  K  are  known  as  the  turbulent  diffusivities  in 
X  y 

the  x  and  y  directions,  respectively.  Then  equation  (4.2)  can 
be  written  as  : 
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Generally,  the  turbulent  diffusivity  is  much  greater  than  the 

molecular  diffusivity,  such  that  can  be  neglected  with 

respect  to  K  or  K  in  equation  (4.4).  Thus,  the  equation  to 
x  y 

be  solved  is  : 
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A  summary  of  the  assumptions  is  as  follows  : 

1.  Two  dimensional  flow  field., 

2.  Negligible  mean  velocity  in  the  verti cal  .di recti  on  ,  v. 

3.  Mean  velocity  in  the  hori zontal ,di recti  on ,  u,  not  changed 
by  the  presence  of  the  gas. 

4.  Gradient  flux  model  for  the  turbulent  dispersion  terms. 

5.  Molecular  diffusivity  negligible. 

The  severity  of  each  of  these  assumptions  will  be 
discussed  in  turn. 


The  first  assumption,  that  of  a  two  dimensional  flow 
field,  infers  that  the  tunnel  can  be  considered  infinite  in 
width.  Thus,  it  is  best  at  the  center  of  the  tunnel,  farthest 
away  from  the  walls.  That  is,  the  walls  induce  three  dimensional 
effects.  In  fact,  the  measurements  were  taken  very  near  the 
center  of  the  tunnel.  Unfortunately,  there  was  no  time  to  make 
measurements  to  verify  the  two  dimensional  character  of  the 
flow. 

In  assuming  that  v  can  be  neglected,  the  term  dropped 
from  equation  (4.1)  is  V3c/3y.  As  3c/3y  can  be  very  large, 
this  assumption  is  actually  not  very  good.  It  is  expected  that 
v  will  become  zero  some  short  distance  from  the  wall,  sc  the 
assumption  causes  large  errors  near  the  wall,  and  that  error  de¬ 
creases  as  the  distance  from  the  wall  increases.  The  actual 
value  of  v  was  not  measured  anywhere  in  the  field. 

The  real  effect  that  the  gas  has  on  the  mean  horizon¬ 
tal  velocity,  u,  can  be  clearly  seen  in  the  experimental  measure¬ 
ments  shown  in  section  3.3.2.  Again,  the  assumption  is  not 
valid  very  near  the  wall,  but  fairly  good  a  small  distance 
from  the  wal 1 . 

There  is  no  indication,  a  priori,  as  to  the  applica¬ 
bility  of  the  gradient  flux  model  for  the  turbulent  dispersion 
terms  in  this  case.  The  assumption  states,  simply,  that  the 
diffusion  terms  can  be  expressed  as  a  function  of  the  gradient 
of  the  average  concentration. 

The  molecular  diffusivity  of  some  typical  pollutants 
is  shown  in  table  4.1.  From  measurements  (see  Ref.  1),  the 
turbulent  diffusivity  in  the  atmosphere  is  usually  on  the  order 
of  1  m2/s.  Hence,  this  assumption  is  quite  good. 


18  - 


4.2  Discretizations 


The  first  step  in  the  discretization  of  equation 
(4.5)  is  the  representation  of  3c/3t  : 


_£  +  l  _£ 

3 c  .  cij  “cij 


(4.6) 


To  use  the  "delta"  representation  (Ref.  3),  we  introduce  the 
parameter  5 . .  : 
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Equation  (4.5)  can  thus  be  written  as  : 
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The  solution  is  sought  for  the  steady  state;  therefore,  a  fully 
implicit  scheme  may  be  used  : 
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Equation  (4.9)  can  be  written  in  operator  notation  : 
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Note  from  equation  (4.7) 
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Substitution  of  equation  (4.11)  into  (4.10)  and  rearrangement 
yields  : 
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The  left-hand  side  of  equation  (4.12)  can  be  factored  after  the 
manner  of  Beam  and  Warming  (Ref.  3)  with  an  error  of  order  At2 
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(Note  that  RHS  refers  to  the  right  hand  side  of  equation  4.12) 
Equation  (4.13)  can  then  be  solved  in  two  separated  parts  : 
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and 


1  +  At 

—  <5-K  i-) 

ax  *  ax 

fij  ■  RHS 


(4.15) 


Equation  (4.15)  is  solved  first,  for  f^;  then  equation  (4.14) 
may  be  solved  for  5^. 

The  diffusive  terms  are  discretized  in  the  following 
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If  h  is  a  stretching  function  that  describes  the  position  y, 
y  =  h (y )  (4.17) 


where  a y  will  define  an  equally  spaced  grid.  By  these 
def i ni ti ons 
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Using  the  same  notation 
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Substitution  of  equations  (4.18),  (4.20)  and  (4.21)  into 
equation  (4.16)  yields  : 
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Each  diffusive  term  is  discretized  in  like  manner. 

The  following  procedure  is  used  to  discretize  the 
convective  terms.  The  stretching  function  is  represented  by 


(4.23 


x  =  g(x) 


and  behaves  in  like  manner  as  h(y).  Thus 
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A  scheme  is  utilized  that  is  partially  centered  and  partially 
excentered,  using  a  proportional ity  factor,  A  : 
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The  convective  term  in  f^j  is  discretized  in  the  same  way. 
Making  the  appropriate  substitutions,  the  fully  discretized 
forms  of  equations  (4.14)  and  (4.15)  are  : 
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4.3  Grid  generation 


The  computational  field  is  a  rectangular  slice  in  the 
two  dimensional  plane.  The  boundaries  at  y= 0  and  y=l  represent 
the  upper  and  lower  walls  of  the  wind  tunnel.  The  boundary  at 
x=0  is  the  inflow  boundary,  taken  a  distance,  b  (which  is  vari¬ 
able),  upstream  of  the  plate  source.  The  boundary  at  x=l  is 
the  outflow  boundary,  taken  a  large  distance  (also  variable) 
downstream  of  the  plate. 

In  the  vertical  (y-axis)  direction,  the  concentration 
gradient  is  most  pronounced  at  the  injection  side,  i.e.,  near 
y=0.  Therefore,  a  finer  mesh  is  required  there  than  in  the  rest 
of  the  field.  Let  the  stretching  function  be  represented  by 
h(y),  where  y  is  the  transformed  coordinate  system.  A  para¬ 
bolic  stretching  function  was  chosen  as  follows  : 


h  =  2  y2/3  +  y/3 


(4.29) 


In  this  way,  the  step  size  at  y=0  is  one-fifth  that  at  the 
lower  boundary  (y = 1 ) . 

In  the  horizontal  (x-axis)  direction,  there  are  sharp 
gradients  at  the  two  plate  edges.  Let  g(x)  represent  the  stretch¬ 
ing  function  in  the  x-direction.  A  cosine  function  was  chosen 
from  x=0  to  the  plate  trailing  edge.  Then  a  parabolic  function 
is  used  from  the  trailing  edge  to  the  end  of  the  field  (x=l). 

The  resulting  function  is  as  follows  : 


b+350 
L 

(4.30) 

b+350 
L 

where 

L  =  the  length  of  the  computational  field  in  mm 

b  =  the  position  of  the  plate  leading  edge  in  mm 

a,d,m,e,f  =  constants  determined  by  the  desired  range  of  step 
sizes  and  total  number  of  steps. 

The  constants  were  chosen  such  that  the  step  size  in  the  center 
of  the  plate  is  five  times  that  at  the  edges.  The  grid  itself 
is  shown,  with  various  parameters  (indicated  on  the  figures)  in 
figures  4.1  to  4.3. 


On  both  upper  and  lower  walls,  the  basic  equations 
for  the  inner  points  are  solved.  Since  centered  discretizations 
in  the  vertical  direction  are  used,  the  basic  equations  use 
non-existent  points.  The  value  of  the  concentration  at  the 
non  existent  points  can  be  replaced  by  the  value  of  the  concen¬ 
tration  at  the  points  in  the  inner  field  that  are  the  "reflec¬ 
tion"  of  the  non-existent  points,  thus  making  a  statement  of 
wall  impermeability.  With  this  assumption,  the  equations  as 
discretized  in  section  4.2  can  be  solved  on  the  boundaries. 

At  the  porous  wall  section,  the  concentration  is  set 
at  a  constant  value,  c  =  0.5.  This  value  was  chosen  as  the 
porous  plate  has  a  porosity  of  50%.  Thus  c  =  0.5  over  the  length 
of  the  plate  should  be  an  accurate  description  of  the  mass  flow 
across  the  plate,  distributed  evenly  across  the  source  plate. 
Thus,  it  should  be  accurate  on  a  global  scale,  but  not  very 
near  to  the  plate  itself. 

4.4.2  Inflow  and  outflow 

For  the  inflow  boundary  (x  =  0)  ,  a  strong  boundary  con¬ 
dition  is  used,  i.e.,  c  =  0.  This  poses  no  problems  if  b,  the 
distance  of  the  plate  leading  edge  from  the  computational 
boundary,  is  taken  large  enough. 

The  boundary  condition  applied  at  the  outflow  (x=l) 
is  on  the  first  derivative.  That  is, 

3c/3x  =  0.  (4.31) 

This  is  only  true  at  a  sufficient  distance  downstream  from  the 
plate  trailing  edge. 


4.5  Test  case 


To  check  the  scheme,  a  problem  is  solved  for  which 
the  exact  solution  is  known.  This  is  done  by  modifying  the 
field  slightly.  The  modified  case  is  one  where  the  plate  ( 
extends  until  x  *  1.  Then  from  equation  (4.5)  : 


3C,  3  3  tt  3  C  ,  3  m  3  C 

—  +  —  UC  =  —  K  —  +  —  K  .  — 


3  t  3  X 


3y  J  3y 


(4.5 


Assume  that  : 


Kx  ■  S  * K 

Steady  state,  i.e.,  3c/3t  *  0 
u  =  1 

K  =  constant,  i.e.,  3K/3x^  =  0 


Then  equation  (4.5)  is  reduced  to 


i£  =  K  +  K 


(4.3 


If  K  is  small  with  respect  to  unity,  then  : 


I  u 


-  27 


Equation  (4.34)  can  be  solved  analytically  with  the  following 
boundary  conditions  : 

x  =  0  implies  that  c  =  0  for  all  y 

y  =  0  implies  that  c  =  1  for  all  x  >  0  (4.35) 

y  =  «  implies  that  c  =  0  for  all  x  >  0 

Note  that  the  last  boundary  condition  can  only  be  used  if  there 
is  no  interaction  between  the  concentration  field  and  the  far 
wall  (at  y  =  1).  The  experimental  data  shown  in  chapter  3 
confirm  that  this  is  an  acceptable  boundary  condition  for  the 
range  of  K  in  which  the  solution  shall  be  made.  Solving  for  c  : 

y/yMx 

c  -  1  -  -2—  [  e  n  dn  (4.36) 

y/ii  > 

0 

or 


(4.37) 


Some  lines  of  constant  concentration  are  shown  in  figure  4.9 
for  K  =  0.01.  The  numerical  program  was  modified  to  use  the 
boundary  conditions  (4.35)  and  the  resulting  solution  is  also 
shown  in  figure  4.4.  It  can  be  seen  that  the  only  discrepancy 
between  the  two  solutions  is  at  the  plate  leading  edge.  The 
small  error  seen  here  is  due  to  the  diffusion  in  the  x-direc- 
tion,  which  is  neglected  in  the  analytic  solution.  It  appears 
upstream  from  the  plate,  as  there  is  no  convection  upstream, 
thus  expression  (4.33)  is  not  valid  in  this  region. 

4.6  Evaluation  of  expected  error  in  scheme 


4.6.1  Upwind1n£ 

To  avoid  the  appearance  of  "wiggles"  in  the  solution 
originating  at  the  boundary,  a  partially  "upwinded"  discretiza¬ 
tion  is  used  for  the  convective  terms.  Assuming,  for  simplicity, 
that  the  first  derivative  of  the  stretching  function  is  equal 
to  unity  (i.e.,  constant  step  size),  equation  4.26  can  be 
written  as  : 
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Expanding  and  simplifying 
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(4.39) 


Thus,  the  process  of  upwinding  introduces  an  effective  diffu¬ 
sion,  with  a  coefficient  of  magnitude  aax/2.  In  any  calcula¬ 
tion  using  this  scheme,  the  magnitude  of  this  term  should  be 
evaluated  with  respect  to  the  physical  diffusivity  of  the 
problem.  Obviously,  A  should  be  kept  as  small  as  possible. 

The  minimum  total  diffusivity  that  the  scheme  can 
handle  is  of  the  order  of  10~3.  Thus,  for  accurate  problem 
solutions  the  physical  diffusivity  should  be  of  order  10-3  or 
greater.  That  is,  the  scheme  is  inapplicable  for  very  small 
physical  diffusivity. 


4.6.2  Error  i_ntr£duced_by  the  ^tretcjii ruj_f unction. 

Consider  the  simplified  equation  with  convection  and 
diffusion  ( i . e . ,  K  =  const.,  u  =  1.),  and  stretching  functions 
in  x  and  y  as  previously  defined.  The  basic  equation,  in 
physical  space,  is,  at  steady  state  : 


=  K  a2<:  +  K  a2c 
ax  ax2  ay2 

This  becomes,  in  the  transformed  space  : 
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(4.41) 


If  a  centered  discretization  is  used  for  the  convective  term. 


ac  „  Ci  +  1 , j~ci -1 , j 
ax  2ax 


(4.42) 


Expanded  in  physical  space. 


and  e 


The  stretching  function  can  also  be  expanded,  in  the  transformed 
space,  as  : 


-  •  AX' 


gi+i  *  v4,V  —  9*  +  H0T 


(4.45) 


’2  M 


g1  _ 1  *  gi - Axg !+  —  g.  +  HOT 


and  equation  (4.44)  becomes 


i£  .  g‘  1£  + 


ix2g'g"  a2c 


(4.46) 


Thus , 


1  ac  _  ac  t  Ax2g"  a2c 
g'  ax  ax  2  ax2 


(4.47) 


It  is  thus  seen  that  due  to  the  stretching  function,  a  diffusion 
term  is  added  with  effective  diffusivity  of  magnitude  Ax2g"/2. 
This  can  be  positive  or  negative,  depending  on  the  sign  of  g". 
With  the  stretching  function  as  shown  in  section  4.3,  and  60 
steps,  the  maximum  value  of  this  effective  diffusivity  is  I0-3. 
So,  again,  the  conclusion  is  reached  that  this  scheme  is  not 
accurate  for  physical  diffusivity  less  than  10-3. 
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4.6.3  Approximation  of_vet oci_ty  due  t o_g r a_v i ty_ef  f  e  c  t  s 

In  the  code  as  written,  the  effect  of  gravity  was  not 

directly  modelled,  but  the  coefficient  of  turbulent  diffusivity 

was  considered  to  include  the  gravity  effects.  This  introduces 

some  error,  obviously,  especially  when  D  and  D  are  taken  as 

x  y 

constant  through  the  field.  That  is,  if  the  is  increased  to 
represent  an  increased  diffusivity  due  to  gravity  effects,  the 
fact  that  it  is  increased  in  both  positive  and  negative  direc¬ 
tions  will  result  in  an  error  for  the  physical  modelling  of  the 
problem.  A  simple  approximation  of  the  velocity  that  could  be 
caused  by  the  density  variation  was  made  to  get  an  idea  of  the 
magnitude  of  this  error. 

Consider  a  two  dimensional  "lump"  of  fluid  of  dimen¬ 
sions  dxdy  and  density  px.  If  the  surrounding  fluid  has  a 
density  pg,  the  potential  energy  contained  in  the  fluid  "lump" 
can  be  represented  by  : 

P  =  ( p  i -ps  )  gdy 

where  g  is  the  acceleration  due  to  gravity.  The  kinetic  energy 
can  be  represented  by 

K  =  Plv2/2 

If  all  the  potential  energy  is  converted  to  kinetic  energy,  and 
there  is  considered  to  be  no  drag  forces  or  viscous  effects 
(thus  a  high  estimation  of  the  velocity),  the  energy  balance 
can  be  expressed  as 

(Pi-Ps)gdy  =  p  i  v2/2 

Solving  for  v, 

,  .1/2 
V  a  j2(p1-ps)gdy/p1j 


The  concentration  gradients  are  the  highest  in  the  y  direction 
near  the  porous  plate,  so  this  equation  was  used  to  calculate 
the  velocity  at  several  positions  near  the  plate.  The  maximum 
value  of  v  was  found  at  the  leading  edge  of  the  plate,  where  v 
had  a  value  of  0.08  m/s.  At  the  trailing  edge  v  was  found  to 
be  0.04  m/s.  This  is  on  the  order  of  1%  of  the  free  stream 
velocity,  which  is  the  same  order  of  magnitude  as  the  minimum 
turbulence  level.  Thus,  we  can  expect  that  modelling  the 
gravity  effects  through  the  diffusivity  will  give  reasonable 
results . 

4 . 7  Presentation  of  results 

4.7.1  Nume r i  ca  1_  resets. 

Several  runs  of  the  numerical  code  were  made  to  see 
the  effect  of  the  pertinent  parameters.  Figure  4.5  is  a  com¬ 
parison  of  two  separate  runs  at  two  stations  downstream  from 
the  plate  source.  All  the  parameters  are  the  same  excepting 
the  diffusivity  in  the  vertical  direction,  Dy.  The  effect  that 
increased  diffusivity  in  the  vertical  direction  has  on  the  concen¬ 
tration  profile  can  be  clearly  seen.  The  maximum  concentration 
(near  the  wall)  is  only  minimally  affected,  but  the  concentra¬ 
tion  is  much  higher  for  the  higher  diffusion  case  as  y  in¬ 
creases.  It  is  proposed  that  "instability"  (defined  in  the 
experimental  section)  would  be  reflected  numerically  by  an  in¬ 
creased  vertical  diffusivity.  In  fact,  the  experimental  data, 
as  evidenced  by  figures  3.4  (j)  and  3.5(h)  reflects  the  same 
trend  as  observed  in  figure  4.5.  In  figure  4.6,  the  effect  of 
varying  uin^-nity  directly  is  illustrated.  The  equations  have 
been  non-dimensi onal i zed  and  the  characteristic  velocity  was 
chosen  as  u^nf1-n^ty.  Thus  the  Uinfinity’  non  dimensi  onal  i  zed , 
is  always  equal  to  1  and  the  effect  of  changing  the  main  stream 
velocity  is  reflected  in  the  value  of  the  non  dimensi onal i zed 
diffusivities,  Dx  =  Dx/UX  and  Dy  =  Dy/UX,  where  Dx  and  Dy  are 
the  non  dimensi onal i sed  diffusivities,  x  is  the  characteristic 


length  scale,  and  U  is  u -j nf -j n -| ty  Thus*  t0  var>y  the  Physical 
main  stream  velocity,  it  is  necessary  to  change  the  values  of 
the  diffusivities  in  both  directions.  This  has  been  done  to 
generate  the  curves  in  figure  4.6.  The  difference  between  the 
two  cases  is  greatest  far  away  from  the  wall,  while  in  the 
experimental  data,  the  effect  of  increasing  the  velocity  was 
felt  the  most  very  near  the  wall.  From  this  it  can  be  concluded 
that  the  change  in  main  stream  velocity  also  has  an  effect  on 
the  diffusivity  that  is  not  modelled  in  this  program,  which  runs 
only  at  constant  diffusivity. 

4.7.2  Compar2son_of  ex£eri mental 
and  numeH  cal^  res^u l_ts^ 

In  figures  4.7  and  4.8,  a  comparison  of  the  experimen¬ 
tal  and  numerical  results  is  made.  In  figure  4.7a  comparisons  are 
made  with  the  experimental  data  when  the  turbulent  intensity  is 
changed.  The  diffusion  coefficients  are  normally  assumed  to  be 
proportional  to  the  turbulent  intensities.  Figure  4.8  supports 
this  assumption,  as  both  diffusion  coefficients  are  doubled  when 
the  turbulent  intensity  is  doubled  and  the  agreement  with  experi¬ 
ment  is  fair.  Figure  4.7b  shows  a  comparison  of  calculation  with 
experiment  for  stable  and  unstable  cases.  For  the  unstable  case, 
the  diffusivity  in  the  vertical  direction  is  taken  as  larger  than 
the  diffusivity  in  the  horizontal  direction.  For  the  stable  case, 
Dy  is  smaller  than  Dx  by  a  factor  of  5,  to  be  consistent  with  the 
hypothesis  that  the  effect  of  gravity  can  be  modelled  by  changes 
in  the  vertical  diffusivity.  To  model  the  effect  of  changing  main 
stream  velocity,  for  the  stable  case,  fairly  good  results  were 
obtained  by  simply  halving  the  diffusivities  in  both  directions 
to  represent  a  doubled  velocity,  as  shown  in  figure  4.8.  In  tradi¬ 
tional  heavy  gas  models  using  the  Gaussian  approach,  it  is  quite 
common  to  use  a^,  which  is  proportional  to  the  diffusivity  in  the 
vertical  direction,  increased  by  a  factor  of  5  from  that  which  has 
been  empirically  determined  for  a  neutral  gas  to  account  for  the 
gravity  effects.  These  calculations  show  that  that  is  actually 
quite  a  good  approximation. 


5.  CONCLUSIONS 


The  comparison  of  the  experimental  data  with  the 
numerical  model  revealed  that  even  with  such  a  simple  model 
for  the  turbulent  diffusion  terms,  the  important  trends  can 
be  predicted.  In  particular,  an  increase  in  the  main  stream 
velocity  is  represented  numerically  by  a  decrease  in  both 
diffusion  coefficients.  It  was  also  verified  that  gravity 
effects  can  be  represented  by  a  change  in  the  vertical  dif¬ 
fusion  coefficient.  For  the  "heavy"  gas  situation,  the  ratio 

of  D  to  D  for  good  agreement  with  the  experimental  data 

y  * 

was  found  to  be  on  the  order  of  5,  which  is  the  commonly-used 
value  for  this  ratio.  And,  finally,  an  increase  in  turbulent 
intensity  of  the  flow  is  well  represented  by  an  increase  in 
the  diffusion  coefficients.  The  relationship  is  directly 
proportional,  as  is  traditionally  assumed. 
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6.  SUGGESTIONS  FOR  FUTURE  STUDY 


6 . 1  Experimental 

Experimental  investigation  of  the  inaccuracies  in  the 
use  of  the  concentration  probe  near  the  wall  should  be  made. 
Also,  the  mean  vertical  velocity  component  should  be  measured. 
The  characteristics  of  the  probe  indicate  that  the  experiment 
should  be  made  with  adjustments  in  the  flow  parameters  such 
that  the  measurement  field  is  extended  further  away  from  the 
wall.  Also,  the  two  dimensionality  of  the  flow  field  should 
be  measured. 

6 . 2  Numeri ca 1 

The  diffusivity  as  used  in  this  code  was  constant 
throughout  the  computational  field.  The  modelization  of  this 
diffusivity  with  respect  to  local  concentration  or  local  tur¬ 
bulent  intensity  should  produce  improved  results.  For  the 
case  where  gravity  effects  are  important,  a  way  to  incorporate 
these  effects  into  the  calculation  of  the  vertical  coefficient 
of  diffusion  should  be  formulated. 

The  numerical  code  could  be  extended  to  three 
dimensions.  This  would  mean  solving  the  factored  equation  in 
three  sweeps  instead  of  two. 

A  similar  program  can  be  used  to  solve  the  transport 
equation  for  second  and  higher  moments,  to  achieve  the  final 
result  of  the  probability  density  function  of  the  concentration 
at  each  point  in  the  field. 
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TABLE  3.1  :  TURBULENCE  MANIPULATOR  CHARACTERISTICS 


FREQUENCY  SHIFT  OF  THE  VERTICAL  BEAMS  =  21  MH 


PINHOLE  DIAMETER  =  50  pm 

DIMENSION  OF  THE  PROBE  VOLUME 

horizontal  beams  Ax  =  .268  mm.  Ay  =  .268  mm, 
vertical  beams  Ax  =  .268  mm.  Ay  =  .268  mm, 

NUMBER  OF  FRINGES  INSIDE  THE  PROBE  VOLUME 
horizontal  beams  25 
vertical  beams  63 

INTERFRINGE  DISTANCE 

horizontal  beams  10.54*10"6  m 
vertical  beams  4.2»10“6 

BAND  PASS  FILTER 

horizontal  velocity  component  30  KH  -  1  MH 
vertical  velocity  component  1  MH  -  3  MH 


Az  =  8.9 
Az  =  3.6 


TABLE  3.2  -  PARAMETERS  OF  THE  LDV  OPTICAL  SYSTEM 


LASER  SOURCE  15  mW  SPECTRA  PHYSICS 
1°  POLARIZER 


CONCENTRATION  PROBE 


VELOCITY 


CONCENTRATION 


FIG.  3. 3- MEASUREMENT  CHAIN  AND  DATA  ACQUISITION 
SYSTEM 
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S  CONCENTRATION  MEASUREMENTS  VS.  HEIGHT  FOR  UNSTABLE  CASE 
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FIG.  3.8-MEAN  VELOCITY  PROFILES  4  =  15. uoo=3<n/s 
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FIG.  4-6 -THE  EFFECT  ON  CONCENTRATION  OF 
A  CHANGE  IN  Ujnf 
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FIG.  4.7  (a)  -  EFFECT  OF  TURBULENT  INTENSITY 
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